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ABSTRACT 

Computerized Tomography (CT) has been demonstrated 
to be a good technique for measuring point-density (void- 
fraction) in two-phase flow systems. Recently, improvements 
have been suggested regarding the choice of filter functions 
m CT methods. These methods are based, essentially, on the 
discrete implementation of the Radon Inversion Formulae which 
is widely used in the medical imaging area. Such methods do 
not require any information, a priori, regarding the distri- 
bution of the density (or the void-fraction). 

A very simple method involving the tomographic chord- 

segment inversion, has been developed and tested for two— phase 

flows having radially symmetric density distributions. This 

method is much simpler and consumes less CPU time relative to 

the more general methods of tomographic reconstruction. For 

test functions the reconstructed density distributions are almost 

exact. For ant air-water bubbly flow data the reconstructed 

3 

values have a maximum deviation of + 0.03 g/cm • The range 

3 

of investigation of the air-water flow data was 0.6 - 0.9 g/cm , 
i.e. void-fraction range of 40% to 10%. These results are 
comparable to the results obtained by the more general methods 
based on the Radon Inversion Formulae. 


1 


CHAPTER-1 


INTRODUCTION 


Computer Aided Tomography (CAT) is being widely 
used in the medical area for the diagnosis of various 
cancerous tissues. The methodology incorporates scanning 
of the patient with gamma rays using appropriate tomogra- 
phic algorithms to reconstruct the density distribution 
of tissues [1] . A basic form of CAT was used in Japan 
in 1946 named "Rotation Radiography" . In this method 
the patient was placed on a rotatographic table, X**»ray tube 
and film were rotated around the patient from 0°-360° while 
the pictures were taken and collected information regar- 
ding various cross sections in the range 0 °- 360 ° . 

This concept of measuring density distribution was 
first investigated by Schlosser et al [2] for a two-phase 
air-water system. The results obtained in void fraction/ 
density measurements have been summarized by Kulacki et al 
[3j. The technique for measuring density distribution has 
great significance because accurate measurements of density 
for various flow systems and components in nuclear systems 
facilitates the computation of heat transfer rates. This 
information is vital from the reactor safety view-point 
because it helps in predicting core burnout etc. This 
technique can be applied in other fields like chemical 
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industries, food-processing and several other research 
areas . 

The various reconstruction methods can be broadly- 
classified in the following categories ; 

(a) Series Expansion Methods 

(b) Transform Methods . 

In the Series Expansion Methods (SEMa) the pixel- 
wise distribution of the function (under-investigation) 
is assumed and then suitable iterative and noniterative 
procedures are applied to achieve the reconstruction of 
the function in the region of interest [4], The iterative 
SEMs are Algebraic Reconstruction Technique (ART), Simul- 
taneous Iterative Reconstruction Technique (SIRT), etc. 

The non-iterative SEMs are Angular Harmonic Decomposition 
(AHD) and Polynomial Decomposition (PD). 

The transform methods are based on the analytic 
formulas based on the Radon Inversion Technique. The 
transform methods are of two types [5] : 

(a) Direct Fourier Inversion (DFI) 

(b) Convolution-Back-Projection (CBP). 

In DFI method the direct Fourier transform of the 
projected data is taken and subsequent 2-d Fourier inversion 
leads to the reconstruction of the unknown distribution. 

In CBP method the data is convolved with a suitable filter 
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function and then back projection of the convolved data 
results in the reconstruction of the unknown distribution * 

An important feature of the tomographic methods is 
that the point -density measurements can be made in a non— 
invasive manner without any prior knowledge of density 
distribution. In non-invasive methods the measurements are 
taken in such a way that the system does not get affected 
which is the case with tomographic methods. 

The currently established reconstruction methods 
are mathematically and computationally complex, so the 
present work is an attempt to develop a simple algorithm 
to measure point densities in radially symmetric flow 
distributions. Such patterns are often encountered in 
gas-liquid flows through pipes. This chord-segment- 
inversion (CSI) algorithm has been demonstrated to be an 
extremely efficient method with further processing resulting 
in the radial density maps. The algorithm has been tested 
against some simulated radially symmetric distributions 
representing bubbly and annular flow . distributions. 

Additionally the CSI method has been applied to 
reconstruct density map for the air-water bubbly flow 
data [2,3,6] . The results appear to be comparable with 
the earlier known more complex tomographic methods [3,7,8], 
The CPU time for the CSI technique is much less than that 
in other general tomographic methods,. 
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CHAPTER 2 

THEORETICAL FORMULATION 

In this chapter a brief explanation of how the 
absorption coefficient is related to the source strength 
and detector reading along a particular chord is given * 
The geometry under consideration is the fan-beam geometry. 
The discrete form of this relationship has been explained 
and the chord-segment matrix and its triangularization 
has been introduced. 

2.1 PRELIMINARIES 

The single-beam radiation attenuation phenomenon 
is represented by 

N = N Q exp [- / c n(r,0) ds ] (1) 

where, 

N = detector reading (count/second) 

N q = source strength (counts/second) 
s = path of radiation (ray) 
c = chord along which s is integrated 
p. = absorption coefficient 
r,0 = cylindrical coordinates. 
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Rewriting Eqn.(1), we have 

d = / c tt( r,0) ds (2) 

where, 

d = In (N q /N) (3) 

In the fan-beam geometry (i.e. the geometry 
in which the beam diverges from the source), for a parti- 
cular chord (ray), Cg, corresponding to the ray making an 
angle © from OS line (See Fig.1), the data is denoted by 

d 9° 

Thus , 

dg = / c ^\i ( r,0) ds (4) 

Here, d^, for different values of-0, is the data to be proce- 
ssed by the tomographic algorithms (in our case the tomo- 
graphic algorithm is CSI). The reconstruction of the 
absorption coefficient, is done and the density,<P;>- , 

(or void-fraction <a>) is determined by a calibration of the 
"CT numbers" (in this case p, ) using some known density 
distributions. The reconstruction has been done for some 

■p •"•I* 

known p, like 1, r, e , e (See Fig. 2). In other words, 
if a set of data dg for any density distribution is known 
to us then we can obtain the CT numbers for that distribution 
and hence by calibration the density. 



j th ring 


( j +1 ) th ring 


— ► S^j - BC-B C 
S - Source 
0 - Object centre 
SO = D 

Fig.1 Data collection geometry 



7 


The discrete form of Eqn.(2) can be written as 

( Fig.1 ) , 


d 


k 


m 

2 S, . p., k = 1,2,. ...m 

j=1 3 


(5) 


•where, 


S 





length of the segment of the kth ray 
falling in the jth ring 
(The hatched lines in Fig.1) 

BC - B ’ C* 

average value of p in the jth ring 
number of rings assumed within the 
object. 


We note that radial symmetry is assumed and 
is now a function of r only. We also note that 


S k,j = °» for 3 < k (6) 

Since the kth ray does not intersect the jth ring if k 
is less than j. Eqn.(5) can be rewritten in matrix 
notation, as 

[d]= [S] [>] (7) 

where, 

[d]= (c^ d m _^ d 1 ) , the data vector, 

M *V_i 


) , the p vector , 
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and 


[S] 


m,m 


m-1 , m 


1 ,m 


0 

S m~1 Jui-1 


5 2,m-1 


0 


0 


is the chord-segment matrix which happens to be lower- 
triangular in this case. 

By Eqn,(7), we get 


trf = [sr 1 [a] 


( 8 ) 


Since the inverse of chord-segment matrix [S] is involved, 
this method is known as chord-segment-inversion technique. 
This results in \i for various rings (or various intervals 
along the radius), A finer data vector will result in 
a better approximation of H along the radial line. 

The expression for S u . ’ s is given by ( See 

K » J 

Appendix~A) 


S k,j = BC ~ B,C1 (See Fig. 1 ) 

= 2rf ( Sin e j+1 ) 5 -( Sin © k ) 2 (Sin ©^-(Sin © k ) 2 


( 9 ) 
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2.2 CHORD-SEGMENT -INVERSION (CSI) METHOD 

The present chord-segment— inversion technique 
is relatively simple and faster than the general tomo- 
graphic methods. 

A FORTRAN program for the CSI algorithm has been 
written and implemented with the flexibility to change 
various geometrical variables, step size, error of inte- 
gration, number of rings (See Appendix I) etc. 

Now, for simulation studies, we want to obtain 
the data vector, [d], from Eqn»(4) . Since we have assumed 
that p.(r,0) is radially symmetric function, i.e. 

1 x (r,0) = p (r) (10) 

do Eqn*(4) can be rewritten as 


d g = ^ ( r) ds . 


( 11 ) 


For simplicity we replace the variable s by x 
and take the origin for x at the mid point of the chord. 
Thus, 

dg = c p. ( r) dx (12 ) 

© 

If for a chord at an angle ©, x.j is the lower limit and 
x 2 is the upper limit for the variable x, then Eqn«(l2) 
can be written as 
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d g = * 2 P- ( r) dx , (13) 

X 1 

By the geometry (Fig„l) and Fig.A1 (Appendix-A) 

— - — . — 2 

r = v (D Sin 9) + x (14) 

Eqn*(l3) reduces to 

x ? 

d D “ / ' p. ( x) dx (15) 

*1 

The steps for the reconstruction of the density 
distribution are as follows s 

(1) Read data [d] in form of a column vector for all 
rays in the fan-beam. 

(2) Compute elements of the lower triangular [S] 
matrix usingEqn,(9) . 

(3) Compute p. -values along the radial segments 
using Eqn a (8) . 

(4) Calibrate Jl^ -values to the density values. 

For the Back substitution, averaging plexi glass 
contribution and data see Appendices B,C,D and E . 

We note that for simulation studies the data vector [dl 
will have to be generated by Eqn*(15). 
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CKAPTER-3 

VALIDATION AGAINST SIMULATED DATA AND RESULTS 
" — — ~ 

3 .1 VALIDATION AGAINST SIMULATED DATA 

In this chapter, we discuss the results for the 
simulated data. We assume the radius of the pipe to be 
one unit and the distance of the source from the centre 
of the pipe to be two units. The pipe is further divided 
into twelve annular rings. 

The CSI algorithm has been tested on the following 
assumed symmetric dis tributions: 

V (r) = 1 .0 

P- (r) = r 

\i ( r) = exp ( r) 

p ( r) = exp ( -r) . 

For above mentioned test functions the errors in recons- 
truction^ 1 0~ 6 , 0,0088, 0,03 and 0,005 respectively. 

Figure 2 shows the results of reconstruction along with 
the actual function distributions. For a listing of output 
see Appendix F , 

The CSI algorithm has also been tested for annular 
flows (See Fig. 3), having the following two distributions: 


QJZZ 


t 


1 



ns 



The reconstruction errors are almost negligible (See 
Appendix G) 0 

The above mentioned-distributions represent ( in a 
calibrated sense) the various types of density/void 
fraction distribution encountered in radially symmetric- 
bubbly and annular flows. 

The reconstruction p, -values matched the assumed- 
values very well for the simulated object of unit radius. 


3o2 RESULTS FOR BUBBLY AIR-WATER FLOWS 

Now here we will discuss the results for bubbly 
air-water flows « The data is taken from the study of 
Ref. [3, 6], Five different data-sets for four different 
cases of density (or void fraction) were processed by the 
CSI algorithm. The algorithm output, ^CTN^, had to be 
calibrated to obtain the density ynhie » For this purpose, 

the previous work, included projection data for a few known 
cases of average density. Figure 4 (and Table FI, Appendix F) 
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shows the calibration chart resulting from the processing 
of the projection data for four known cases of densities. 

The output, CTN, was corrected to the base of air, CTN = 0.0, 

t 

to eliminate the effect of plexiglass* (Details in Appendix D). 

The four points are joined by a straight line such 
that the line represents the best fitting line for these 
four points. This figure is now the calibration curve 
because now by knowing the < CTN> we can obtain the 
corresponding density value . 

Similarly the data for all scans (See Appendix £.) 
has been processed and we get the output (See Appendix H) 
in such a format that it gives the value of CTN for corres- 
ponding anglc/radius • The value of P (obtained after calibra- 
tion) for each radius has been plotted to show the recons- 
tructed density and hence reconstructed profiles for 
various void fractions (See Figs 5-9). 

In Figure 10 the comparison of reconstructed 

densities with actual densities has been shown by an 
alternative method [3.6] in which X-axis is the actual density 
and Y-axis is reconstructed density,, A line such that 
<P>- < p qj > has been drawn and various densities obtained 
from calibration method are located to show the deviation 
of the reconstructed densities from the actual densities. 
Appendix F summarises the <CTN> and <P>results for all 
scans in a tabular form. 


gm/cm 







<0> -1-<0d'> 

Legend 

3 

<p>- 0*9 gm/cm 
<P>-0*8 gm/cm 3 
<9> - 0*7 gm/cm 3 
<9> -0*6 gm/cm 3 


r (cm) 

Density profile (radial) for scan 2 




< 9 > = 1 - <<*> 
Legend 

o <£> =0-9 gm/cm 
d <£> -0-8 gm/cm 
A <^> -0*7 gm/cm' 
• <9> =0*6 gm/cm' 


Density profile (radial) for scan 3 





< 9 > 

Legend 

3 

o <p> =0*9 gm/cm 

3 

□ <p> = 0*8 gm/cm 

A <p> -0-7 gm/cnr? 
• <P> =0*6 gm/cm" 


Density prof He ( radial ) for scan 4 




<p>=l-<a> 

Legend 

4 

<p> - 0*9 gm /cm 

<p> = 0*8 gm /cm 

* 

<p> -0*7 gm /cm' 
<p> r 0*6 gm/errf 


r (cm) 


Density profile (radial) for sc an 5 
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The maximum error obtained in reconstruction 
3 

is + 0.03 gm/cm . In Fig. 11 the variation of error with 
density has been shown. 

The maximum relative errors are + 4% and -5% 
respectively. In Fig. 12 the variation of relative error 
with density has been shown. 

We note that there are some statistical fluctuations 
in the count-rate recorded by the detector. This uncertain- 
nity leads to an erroneous reconstruction. Such a discre- 
pancy appears to be quite obvious for constant density 
cases of air, pine, walnut and water (See Fig. 13 ) . 

However assuming Poisson distribution and applying + 1 a 
and ± 3 a corrections (where 1 a implies one standard 
deviation of the count rate, N), the ripple appearing in 
Figure 13 is smoothned out as is evident from Figures 14 
and 15, This exercise leads to a P CT "band" for the air- 
water flow data. Since the point density values were not 
available to make any meaningful comparison , the air— water 
cases have not been presented. 


Error in density measurement for all 
scans 




0-4 


0*6 


1*0 


Fig. 12 


<P> 

Relative error 


0*8 . 
,g/cm' 

for all 


scans 








Fig.13 Calibration cases with mean value 
of counts 




POISSON STATISTICS 
N : Mean value of counts 


o Water <p> = 1-0 g/cm' 5 

A Walnut <p^ =0-73 ” 

x Pine <p> = 0-41 » 

o Air <p> =0-0 


n 1 = n-Jn 

n 2 -N+jN 





— a- - 


'S — — -8~ — -P 




=& 








=■=*=“ 


X 




,\ 




CTN 
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o Water 
A Walnut 
x Pine 
o Air 


<p> - 1-0 glcrr? 
<f> - 0-73 v 
< p> - 0-41 
< - 0*0 ,j 


POISSON STATISTICS 
N : Mean value of count 

--- N 1 = N - 3jN 
N 2 = N+3,/R 
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CHAPTER-4 

DISCUSSION 


The proposed chord-segment inversion (CSI) technique 
has been demonstrated to be applicable to radially symme- 
tric two-phase flows. The results are comparable to the 
results obtained by more general methods of tomographic 
reconstruction [2, 3,6, 7, 8], 

However, taking the advantage of the radial symmetry 
in the theoretical formulation itself, the CSI method 
consumes less Central Processing Unit (CPU) time and it is 
simple to use. 

Comparable results have been obtained by CSI 
technique for density measurement in bubbly air— water flow 
for cases 0.6 <_ p < 1,0, The maximum deviation of P q j 
is about + 0,©3g/cm which is approximately +5%. The 
maximum relative errors are +4% and -5% respectively. 

The error in the point density values could not be estimated 
as no alternate method was available in the study of Ref.[3,6], 
But the test-function results indicate that for ’’pure" 
data the reconstruction is exact. The results of CSI method 
are also comparable to "radial polynomial" method already 
suggested and tested [9] for radially symmetric flow 
distributions. 
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APPEND IX-A 


Suppose we have a pipe of radius R and a gamma 

ray source S at a distance D from the centre of the 

pipe. Consider now jth and (j+l)th annular rings and 

five chords , C Q , C Q , C Q and C Q such that 

y o y k y j y j+1 y max 

they make angles 0 3 0,. 0 . , © . , , , © as shown in 
1 J k 7 j J j+1 max 

Fig.1 , 


The angle corresponding to chord is the 

y max 

maximum angle © and is given by 


Sin(0 max ) = R/D (A ~ 1) 

©max = Sin " 1 (R/D) (A “ 2) 

1.U 

Note that any ray (chord) C Q is the tangent to the (n-l) x ' 

y n 

annular ring. 

Now we want to calculate the length of the segment 

of the kth ray falling in the jth annular ring that is S, 

K »3 

From Fig, 1 this length is the hatched line 

S k,j = BC - B ' C ' (A~3) 


Fig.A1 , 


Now see the detailed geometry of Fig, 1 in 
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Figure A1 : Fan Beam Geometry Details 


We are interested in calculating the length LL’ in the 

Cq chord such that C Q is the tangential chord, 
k to max 

From the triangle SOL’ 


D R 

SirKin Sin(© k ) 


( A— 4 ) 


p = sin -1 (2_S|nLSid 


But from Eqn.(A-l) 


D Sin (© mav ) 
max 


So from Eqn.(A-5) 

P 


1 D Sin( © v ) 

sin ( D - sinCeJ" ) ) 


max J 


(A-5) 


(A— 6) 


(A— 7 ) 
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, Sin (© k ) 

Sln ( 's'i"re~J ) 


(A-8) 


The length LL* 


LL* = 2R Cos(|3 ) 


( A-9 ) 


_< Sin(0 v ) 

= 2D Sin(g ax ) Cos(Sm \ ?r j; r £ r ^ ) ) (A-10) 


max - 


Sin 5 © 


20 Sir *W 1 - TTS 


Sin ©. 


max 


(A— 1 1 ) 


2D * Sin © ~ Sin^9 k 


max 


( A— 1 2) 


By this token we have in Fig* 1 




r“2: 


BC = 2Df Sin 0 . +1 - Sin © k 


( A— 1 3 ) 


. — - 2 

B»C f = 2D^ Sin ©j - Sin © k 


(A— 14) 


hence the length 


s k,j = BC ~ B,c * 


= 2D[' / ‘sin 2 ©j +1 -Sin 2 © k Sin 2 ©^ - Sin 2 © k ) 1 


(A-15) 
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The Back' Substitution. is as follows: 
Recalling equation (7) 

[ d]= [S] [uJ 

The expanded form of this equation is 


(B-1) 





1 1 

— 


S 

m,m 


0 ! 

V m 
m 

■ 

S , 
m-1 ,m 

S j. S A 

m-1 ? m— 1 

iT 


1 

m-1 

V 

O 

0 

O 


« 

0 

0 



• 


Q 

© 

0 



J 

S 4 

1 jin 

S 2.m-1 

S 1,1 


“1 


_ 


1 



(B-2) 


Then mho algebraic equations will be 


m 


^m-l 


'm-2 


S ’’ 
m.m -m 


5 m-*1 ,m ; S m-1 ,m~1 ^ m-1 


( B— 3) 
( B— 4) 


^m~2,rahm *" ^m-2*ni-*1 ^ m-1 + ^m~2 , m-2 ^m-2 

0 0 0 ( B— 5) 


Si U- + S,-. jiP' *+ 9 » a t S> a P" 4 

1,m r m 2,m~1 m-1 0 *° 1 


(B-6) 


From Equation (B-3) we have, 


Substituting this value of yt in (B—4) we obtain the 
value of and then using the value of in 

next equation to get ^ m _ 2 and so on. 

So by back substitution we get the values of all 



APPENDIX C 


t ot the simulated data study of radially symmetric 

distributions, the p has been reconstructed. In Fig. Cl, 

Pj an< 3 are the reconstructed values of p ’s 

corresponding to radius 0, R and R respectively. Simi- 

J «J * ^ 

larly p ^ , p ^ and p j+1 ore the actual p’s at radius 
0, R| and respectively. . '7 JJj + ^ 



R q ,Cl AVERAGE- VALUES _ _ 

^ Now take the average values of p’s < p , <p > 

J * 3 

and < p such that they are the values of p’s in first, 


second and third annular rings 


Now < p ^ > is the value 


of p at radius R Q as shown in Fig. Cl 



If area between R. +J) and R q is tnen 


A 


1 




Similarly if area between R and R. is A 0 then 


A 

Taking A 


2 

1 


= it R 




R 


o 




+ 



i.e. the average radius of the annular ring having R. 
internal radius and R^ + ^ outer radius. This must be taken 
into account when calculating |i at various radii. 


APPENDIX D 


The data for all scans has been taken from [6] • 

There was plexi-glass around the pipe in the experiment 
of [6] for all scans. For air the CTN (absorption 
coefficient) must be zero but for air data set we do not 
get the CTNs ( corresponding to different angles) equal 
to zero but some non-zero values, this is due to the plexi- 
glass whatever CTNs we get will be the CTNs for plexi- 
glass but not for air. 

So this plexi-glass contribution must be taken into 
account in order to obtain the correct values of CTNs for all 
cases. For this we make air as the reference and all CTNs 
equal to zero by subtracting CTNs of plexi-glass. Similarly 
for any case we subtract the CTNs of plexi-glass to get 
CTNs and hence correct <CTNs> . 


the correct 



APPENDIX F 


TABLE -FI 


(Results for air-pine walnut and water) 


Case 

<P>, g/cm 3 

<CTN> 

Air 

0.0 

0.000 

Pine 

0.41 

0.075 

Walnut 

0.732 

0.137 

Water 

1 .00 

0.201 


TABLE-F2 

(Reconstructed Densities) 




TABLE -F 3 

(Reconstructed <CTNs > ) 




Table F 1 Scan No. 

































































































































































































































Table E.2 Scan No. 
























































































































I able I .3 Scan No. 
































































































Table E.4 Scan No. 






















































































































Table E.6 Scan No. 
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